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Abstract. We study the effect of free boundaries in finite magnetic systems 
of cubic shape on the field and temperature dependence of the magnetisation 
within the isotropic model of D-component spin vectors in the limit D — > 
oo. This model is described by a closed system of equations and captures 
the Goldstone-mode effects such as global rotation of the magnetic moment 
and spin-wave fluctuations. We have obtained an exact relation between the 
intrinsic (short-range) magnetisation M = M{H,T) of the system and the 
induced magnetisation m — m{H,T) which is induced by the fleld. We have 
shown, analytically at low temperatures and fields and numerically in a wide 
range of these parameters, that boundary effects leading to the decrease of M 
with respect to the bulk value are stronger than the finite-size effects rendering 
a positive contribution to M. The inhomogeneities of the magnetisation caused 
by the boundaries are long ranged and extend far into the depth of the system. 
PACS: 75.50.Tt - 75.30. Pd - 75.10.Hk 



1 Introduction 

Small magnetic particles have been of much interest owing to their technological applications, mainly in the 
area of information storage. From the experimental and theoretical point of view, these systems are very 
interesting for they show superparamagnetism at high temperature and exponentially slow relaxation rates 
at low temperature due to anisotropy barriers. Although in most of theoretical approaches to the dynamics 
of a small magnetic particle the latter is considered as a single magnetic moment, deviations from this simple 
picture become crucial with the reduction of the system size. In magnetic nanoparticles (see, for a review, 
Ref. [Ql), the contribution of the surface to the thermodynamic properties becomes comparable with the bulk 
contribution, and the magnetisation and other characteristics may be spatially inhomogeneous. The latter 
is realised due to additional thermal disordering near the surface at elevated temperatures or due to the 
breaking of symmetry of the crystal field for surface spins, which may result in a strong surface anisotropy. 
Another manifestation of the symmetry breaking at the surface is the possible unquenching of the orbital 
moments of surface spins, which may be responsible for a significant increase of the particle's magnetic 
moment, e.g. in 3d elements [Q. 

A finite-size magnetic system with realistic free boundary conditions presents a spatially inhomogeneous 
many-body problem. A great deal of work up to date has been based on the Monte Carlo (MC) technique 
for the Ising model. MC technique was also used with the more adequate classical Heisenberg model to 
simulate idealised isotropic nanoparticles with simple cubic (sc) structure and spherical shape in Ref. 
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Magnetic nanoparticles with realistic lattice structure were simulated in Ref. Q taking into account surface 
anisotropy and dipole-dipole interaction. 

Familiar analytical methods such as the mean-field approximation (MFA) and spin-wave theory (SWT) 
are (at least in their standard form) inapproriate for finite magnetic systems because of the Goldstone 
mode corresponding to the global rotation of the magnetic moment in zero field. A more sophysticated way 
to analytically treat finite magnets is to take the limit D oo for the model of /^-component classical 
spin vectors, which was introduced by Stanley Stanley has shown Q that for Z) — s- oo the partition 
function for this model in the bulk coincides with that of the exactly solvable spherical model (SM) [Q. 
Both models provide a reasonably good approximation (about 5% overall accuracy) to the isotropic classical 
Heisenberg model in three dimensions. For spatially inhomogeneous systems such as magnetic particles or 
films with free boundary conditions, these two models become nonequivalent, the SM having been shown 
to yield unphysical results because of the inappropriate global spin constraint Q . This drawback was fixed 
in an improved version of the SM using local spin constraints on each lattice site ^, |l^ . Yet two profound 
differences between the Z? — > oo model and the SM cannot be eliminated. First, there are two (longitudinal 
and transverse) correlation functions in the former [ll| and only one CF in the latter. Second, generalization 
for the anisotropic case is only possible for the — > oo model. The resulting anisotropic spherical model 
(ASM) was applied to describe phase transitions in domain walls ||l^ and in magnetic films |l^. The main 
feature of the D ^ oo model or the ASM is that it takes into account Goldstone modes in the system, which 
prevent phase transitions for isotropic systems in dimensions d < 2. In contrast with the linear SWT, this 
model self-consistently generates a gap in the correlation functions which avoids the infrared divergencies. 
Anisotropy plays a crucial role in these systems since it breaks the symmetry of the system and makes phase 
transitions in two dimensions possible. 

So far, the ASM was only applied to spatially inhomogeneous systems in the plane geometry |l^, p^ . 
Here we extend it for finite box-shaped magnetic systems with free and periodic boundary conditions (fbc and 
pbc), restricting ourselves to the idealised case of isotropic coupling and simple cubic lattice and reserving 
consideration of realistic lattice structures, bulk and surface anisotropics, etc., for the subsequent work. 
We obtain analytical asymptotes at low temperatures and fields and a complete numerical solution in the 
whole range of T and H for the intrinsic magnetisation M = M{H,T) and the induced magnetisation 
m = m{H,T), and we prove the exact relation m = MB(MMH /T) between them, N = N^NyN^ being 
the number of sites in the system and B{x) the Langevin function. For the system with fbc, the value of 
M is significantly reduced with respect to that for the pbc model and to the bulk due to boundary effects. 
We consider temperature dependences of local magnetisation Mi and estimate the critical indices for the 
magnetisation at the faces, edges, and corners of the system. We show that because of the Goldstone modes, 
the inhomogeneities of Mi extend far into the depth of the particle, i.e., the magnetisation profile in our 
isotropic model is long ranged. This feature is in accord with the MC simulations for the Heisenberg model 
H, U] , as opposed to the MFA predictions j|] . 

The rest of the paper is organized as follows. In Sec. ^ we give the equations describing the spatially 
inhomogeneous ASM in a magnetic field, generalizing the results of Ref. [l^ . Further, restricting ourselves 
to the isotropic model, we define the two magnetisations mentioned above and derive relations between them. 
In Sec. Ijwe consider the model with periodic boundary conditions and derive low-temperature and low-field 
analytical expressions for the intrinsic magnetisation. In Sec. ^ we take into account boundary effects for 
the fbc model at low temperatures. In Sec. || the numerical results for the temperature, field, and spatial 
dependencies of the magnetisation of finite magnetic systems are presented. 

2 Hamiltonian and basic equations 

We start with the Hamiltonian of the uniaxially anisotropic classical _D-component vector model, which can 
be written in the form 



where is the normalized D-component vector, |si| — 1, and 77 < 1 is the dimensionless anisotropy factor; 
H is the magnetic field, and Jy the exchange coupling. 




(2.1) 
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In the mean-field approximation the Curie temperature of this model is T^^^ — Jq/D, where Jq is 
the zero Fourier component of Jij. It is convenient to use T^^-^ as the energy scale and introduce the 
dimensionless variables 

= T/T^^^, h = H/Jo, \, = J,,/Jo. (2.2) 

For the nearest-neighbour (nn) interaction Jij with z neighbours, Ay is equal to 1/z if sites i and j are 
nearest neighbours and zero otherwise. 

Using the diagram technique for classical spin systems |l^, |l^ in the limit D — > oo and generalising 
the results of Ref. for spatially inhomogeneous systems to include the magnetic field h =hr^ex + h^e^, 
one arrives at the closed system of equations for the average magnetisation components a = 1 (z) and a = 2 
{x) and correlation functions for the remaining spin components labelled by a > 3, 

Here all correlation functions are equal to each other by symmetry. The system of equations describing the 
anisotropic spherical model consists of equations for the magnetisation components 

rrixi = Gi{ha; +ri'^Xijmxj), m^i ^ dih^ +^ Xijin^j), (2.4) 
j j 

and the Dyson equation for the correlation function 

Sii = OGiSii -f 'qGi ^ \jSju (2-5) 

3 

where Su is the Kronecker symbol. Gi is the so-called gap parameter to be determined from the set of 
constraint equations on all sites i — 1 , . . . , A/" of the lattice 

s,, + ml = 1, (2.6) 

where mf = m^j -I- m^j. 

The system of equations above describing the anisotropic spherical model (ASM) will be used in its 
full form elsewhere. In the present work, we will consider isotropic (ry = 1) magnetic systems of box shape 
of volume N — N^NyN^, with simple-cubic lattice structure, and nearest-neighbour exchange coupling, 
in a uniform magnetic field. In this model, the magnetisation m is directed along the field h, thus one 
can set h =hez and m,=m,e2. In the limit of Z? — > cxd it becomes immaterial what is the exact number 



of "transverse" components: Z) — 2 for the anisotropic model, Eq. (2.3), Z) — 1 for the isotropic model 
in field (the case we are considering here), or D for the isotropic model in zero field. In general, there 
are a few components with nonvanishing averages (sq.), and the remaining components called "transverse" 
components, are booked into the correlation functions. 

Solving the system of equations above consists in determining and as functions of Gi from the 



linear equations (2.4) and (2.5), respectively, and inserting these solutions in the constraint equation (p.6|) in 



order to obtain Gi. The two main types of boundary conditions for our problem are free boundary conditions 



(fbc) and periodic boundary conditions (pbc). In the case of fbc, if the summation index j in Eqs. and 



(2.5) runs out of the lattice, the corresponding terms are omitted. In this case, m.i and Gi are inhomogeneous 
and Sij nontrivially depends on both indices due to boundary effects. In the pbc case the solution becomes 
homogeneous and greatly simplifies. Although the model with pbc is unphysical, it allows for an analytical 
treatment at all temperatures and for the study of finite-size effects separately from boundary effects. 



One can also introduce a matrix formalism and rewrite Eqs. (2^) and ( p.5| ) for the isotropic model in 
field as follows 

^V^jnij = h, ^'DijSji = dSii, (2.7) 
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where we have defined the Dyson matrix T>ij = G^ ^5ij — Xij. The solutions of these linear equations are 



^hJ2'Dij\ s.j^eVi^, (2.8) 
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T> ^ being the inverse of the Dyson matrix T>. Substituting these solutions into the constraint equation (2.6) 
results in a closed system of nonlinear equations for the gap parameter Gi 



(2.9) 



The average magnetisation per site defined by 



m 



--y 



(2.10) 



vanishes for finite-size systems in the absence of magnetic field due to the Golstone mode corresponding to 
the global rotation of the magnetisation. On the other hand, it is clear that at temperatures <C 1 the spins 
in the system are aligned with respect to each other and there should exist an intrinsic magnetisation. The 
latter is usually defined for finite-size systems as 



M = 



\ 




\ 



1 



(2.11) 



where the second expression is valid in the limit D oo. Here m and Si 



are defined by Eqs. ( ^.3| ) and 
In this case in the limit 6^0 
one has Sij ~ 1 for all i and j, and M —* \. For 6 ^ oo the spins becomes uncorrelated, Sij = Jy, and 



(2.1C). One can see that M > m. Note that M remains non zero for h = 0. 



M 1/^/77. In the limit of A/" — > oo, the intrinsic magnetisation M approaches that of the bulk system. 

The applied field suppresses the global-rotation Goldstone mode and thereby renders the magnetisation 
m of Eq.(2.10) non zero. Therefore, it is convenient to call the latter the induced magnetisation, in contrast 
with the intrinsic magnetisation AI. If the field is strong or the temperature is low, the spins align along the 
field, the transverse correlation functions Sy in Eq. (2.11) become small, and the magnitude of the intrinsic 
magnetisation approaches the induced magnetisation. 

One can establish an important relation between the intrinsic magnetisation M and induced magnetisa- 
tion m. To this end, one first obtains from Eqs. (C 



and (2.10) the useful relation 



1 X ^ 9 m 



M h 



(2.12) 



Then, substituting the latter in Eq. (2.11) and solving the equation thus obtained for m, m? + Om/{Mh) — 
M2 = 0, yields 

m = M "^^Mh/e ^ MB{J7MH/T), (2.13) 

1 + ^1 + {2NMh/dY 



where 5(0 = (2^/1))/ 1 + + 



is the Langevin function for Z? 3> 1. One can find in the litera- 



ture formulae of the type m — MgBiMMsH /T), where Ms is usually associated with the bulk magnetisation 



7T) 



at a given temperature (see, e.g., Refs. |18|, |19|). In o ur ca se, Eq. ( 2.13 ) is exact and M — M{T,H) is ex- 
plicitly defined by Eq. 
h ^ hy where 



( 2.11 ). For large sizes Af, Eq. ( 2.13 ) describes two distinct field ranges separated at 



Moie) = M{0,h = 0). 



(2.14) 



In the range h<hy, the total magnetic moment of the system is disoriented by thermal fluctuations, and 
the induced magnetisation m is lower than the intrinsic magnetisation M. In the range h>hv, the total 
magnetic moment is oriented by the field, m approaches M, and both the latter further increase with field 
towards saturation (m ^ M = 1) due to the suppression of spin waves in the system. This scenario is quite 
general and inherent to all 0{D) models, as was shown on phenomenological grounds in Ref. Early 
Monte Carlo simulations of Ref. for isotropic Heisenberg systems of spherical shape confirm Eq. (2.13) 
within statistical errors, although the accuracy is not high enough to decide whether this relation is exact. 
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One can also consider local magnetisations. The local induced magnetisation is simply the vector rxii, 
while the local intrinsic magnetisation can be defined as follows 



1 



^^^^Mi^'-jfT^^^/^miT. (''^^ + • "^^) (2-15) 

One can check the identity (l/JV) Mi = M showing the self-consistency of the definition given above. 
We start by considering the pbc model in the next section. 



3 Systems with periodic boundary conditions 

In this case, as was said above, the system becomes homogeneous, that is Gi = G, = m, and the 
correlation functions can be found by performing a Fourier transformation of the type 

F. = ^^e-kr.^^^ ^^Y^^^F^^ (3-16) 
k i 

where for a cube (A^ = N"^) one has 

ka = 2i:na/N, ria = 0, 1, . . . , — 1, a — x,y,z. (3-17) 

This results in 

«« = }^E*k = ^GPjv(G), (3.18) 
k 

where 

and Aj^ = Jj^/Jq. In the lattice Green function P/v(G) the contributions of the would-be Goldstone mode 
with k = and other modes have been separated from each other (hence the prime on the sum in Pn{G)). 
In the bulk limit ^ oo for the sc lattice 

P.,C, . P,G, . { , (3,20) 

where d = 3, the Watson integral W = 1.51639, and cq = (2/7r)(3/2)'^/^. The difference between the sum 
and the integral, 

k 

describes the finite-size effect and for the sc lattice, Aj^ = (cos/cj, -f cosfcj, + cos kz)/d, d = 3, behaves as 1/A^ 

Ar'.i^^--5^, (3,22) 

The easiest way to obtain the coefficient in this formula is to plot A^'''^-' vs 1/A^ (see Fig. p. For finite N, 
P/v(G) becomes regular at G = 1. In the large- TV limit for 1 — G ^ l^min ^ V^^; where /cmin ^ 1/iV is the 
minimal wave vector in a finite-size system, one has 

P^(G)-<P'^^)-c^A^(l~G), c^ = |^ E' 0-382. (3.23) 
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For larger 1 — G, the square-root singularity in Eq. ( ^.20 ) is restored. More precisely, one obtains, 



(pbc 
N 



cn N{1~ G), 1 - G < l/N^, 
cnVT~G, 1 - G > 1/7V2. 



(3.24) 



Using \ij = 1 in the second of Eqs. ( p.4[ ), one arrives in the isotropic case rj — 1, at the system of 
equations 

hC 

m= -, n? + eGPN{G)^l. (3.25) 

1 — G 

The solution G{9) of these equations decreases with increasing temperature, and at high temperatures, the 
leading asymptote is G = 1/0. For the bulk system at low temperatures or high fields, G tends to 1/(1 + h). 
For finite-size systems, even in zero field, the Goldstone-mode contribution 1/[A/'(1 — G)] in Eq. ( 3.19| ) makes 
G smaller than unity, 1 — G^9/Ma,t6^ 1, which means that there is a gap in the correlation function 
sq. As a consequence, the magnetisation m in Eq. ( p. 25 ) vanishes for ft. ^ for any finite Af. The same 
happens for low-dimensional systems d < 2, even in the bulk limit. In contrast, for three-dimensional bulk 
systems in zero field, G remains equal to 1 and the gap in the correlation function closes for 9 < 9c ^ 1/W. 
This is the reason to call G the gap parameter. Below 9c the spontaneous bulk magnetisation is given by 



(3.26) 



The intrinsic magnetisation of Eqs. ( ^.11 ) can with the help of Eqs. ( ^.12 ) and the first of Eqs. ( |3.19 ) 
be rewritten as 



M 



9G 



Nil - G) 



1-9GPn{G), 



(3.27) 



where the last equality was inferred from the second of Eqs. ( ^.25 ). Note that m and M are related by Eq. 
(2.13) which is valid for all types of boundary conditions. Under a very small magnetic field, we can write 



M 



B{AfMHlT) 



NMH H 



T 



where Hy = is the field at which the global rotation of the particle's magnetic moment is suppressed. 
The corresponding reduced field hy = y- was defined in Eq. (|1|). Here 1 - G can be found from the first 
of Eqs. (E2|) and Eq. (|2l3|) as 



1 - G ^ h/m = {h,/2)[l + ^l + {2h/Ky 



(3.28) 



With the help of these equations, the dependence M = M(h, 9) can be established perturbatively in different 
regions of the parameters. In particular, below 9c for Af ^ 1 and h <^ h^, one can use Eq. (3.24) with 
1 — G ^ hy[l + [h/hyY], which results in 



2nib 



(pbc; 

N 



N^CN ,2 



h <^ hy. 



(3.29) 



The last term in this expression cannot be obtained from the naive spin-wave theory. 

On the other hand, in high fields, there is a crossover to the bulk spin-wave singularity at the field 

2 T 

Hs = ^jfrjj which is much larger than Hy- Indeed, for hg with 



Jo 



2d 



(3.30) 



Eqs. (3.27), together with the second line of Eq. (3.24), leads to the following high-field behaviour of the 
magnetisation 



M = mb + 



2mb 



(pbc) 
N 



coVhj 



(3.31) 
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which shows the well known singular \fh spin-wave correction to the magnetisation in three dimensions. 

For hy h <^ /is, Eq. (3.27) and first line of Eq. ( |3.24 ) yield the linear field behaviour of the 
magnetisation as a function 



2mb 



-WA 



(pbc) 
N 



+ NcNh 



(3.32) 



Therefore, we obtain two crossovers in field, one at hy given by Eq. (2.14) between the quadratic and 
linear behaviour, and another one at hs defined in Eq. ( 3.30| ) between the linear and square-root behaviour 
for the magnetisation. 

Note that both of these crossovers occur for the induced magnetisation to, Eq. ( 2.13| ), as well as the 
intrinsic magnetisation A/, Eq. (3.27). At low temperatures in zero field M deviates from 1 according to the 
law 

M = 1 - 9wj^^"'''/2, (3.33) 



where the coefficient in the linear-6' term is smaller than in the bulk [see Eq. ( 3.22| )]. 

It should be noted that in his early work on thin magnetic films Doring [^0[ used expressions of the type 
of Eq. ( 3.33| ), where the dangerous k = mode was excluded on physical grounds. This intuitive approach 
was applied to clusters of quantum spins in Ref. pl| using the numerical solution for the energy levels in 
the linear spin-wave approximation. Whereas the results of these works for the temperature dependence 
of the magnetisation are reaso nable if one specifies the meaning of the "magnetisation" as the intrinsic 
magnetisation in zero field, Eq. ( 2.13 ) relating the intrinsic and induced magnetisations with each other and 
describing global fluctuations cannot be obtained if one just neglects the k = mode. 

In spatial dimensions d < 2, P{G) diverges for G ^ 1, which rules out the long-range order in bulk 
systems. For any finite-size system, however, M — > 1 at low temperatures. In particular, in two dimensions 



one has W 



(pbc 

N 



Stt In + const, which results in 



A/ = 1 - 6'(47r In + const), A^ > 1, 



< 1. 



(3.34) 



F or th e fields h ^ 1,6 one has, again, to ^ 1 , G = 1/(1 + /i) <C 1, P{G) ^ 1, and thus using the first of 
Eqs. ( 3.25 ) and the second equality in Eq. ( ^3.27 ) one obtains 



M = 1 - 



2{l + h)' 



(3.35) 



Here keeping 1 in the denominator improves the asymptotic behavior for moderately high values of h. 

For the model with free boundary conditions, one cannot obtain a general analytical solution, and in the 
main range of temperatures and fields the problem has to be solved numerically. Analitycal solutions exist, 
however, for high and low temperatures and for high fields. Before proceeding with the numerical solution 
of our problem, we will consider in the next section the most interesting analytical solution of the fbc model 
for 6* < 1. 



4 Boundary effects at low temperature 

At zero field the induced magnetisation m is zero, but in the limit of zero temperature all spins in the system 
beco me s trongly correlated with each other: Sij 1. For <C 1 one can search for the solution of Eqs. (2.5) 
and (^.6|) in the form 

(4.36) 



I -Si 



G,; = Gf ^ - dG,, 



where Gj-"^ is the zero-temperature value of Gi, Ssij and SGi are small corrections. The zeroth order of Eq. 
( |2.5| ) becomes 

l = Gi'^Y.^^^ (4-37) 

3 
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which determines G^f^ . If the site i is not near the boundary so that all sites j are within the system, then 
the sum above and thus Gj-"^ are equal to one. If i is near the boundary, some of the values of j run out of 
the system and are not counted, thus G^^ increases. For the sc lattice one has 

^ =^-l{(>^^,l+S^,,N+S^^,l+Si^,N + Si^,l+S^^,N), (4.38) 



(-(0) 6' 

where i = {ix, iy, iz}- That is, G-°^ = 6/5 on the faces, G-°'' = 3/2 on the edges, and G^°^ = 2 at the corners. 
To first order in 6, Eqs. (p!q) and (^.7[) become 



Ssii = 0, 

J2v[fSs,i = -eSu + 6GJ{G^f^)\ (4.39) 

3 

where Pj-^"* is the Dyson matrix defined in Eq. (2.7) taken at zero temperature. The solution of the linear 
problem above has the form 

<5G. = ^(Gfy5:p(°)p-(°) 

i 

Ss,,^e[V;,'^"^-V;^'^% (4.40) 

where V^j^^^"^ are the matrix elements of the inverse Dyson matrix I?'^*^-'. 

The latter equations can be used to calculate different observables. In particular, the intrinsic magneti- 
sation M of Eq. (2.11) becomes 

To evaluate the coefficient in the linear-6' term, it is convenient to expand T>^^^^'^ over the set of eigenfunctions 
. of the Dyson matrix vf^ : 

^-1(0) ^^:^W^^ (4.42) 
k 

where the eigenfunctions satisfy 

E PkM? = Mk^^k, ' E ^k.^k'. = -^kk' • (4-43) 

i i 

The latter can be found exactly since the matrix P^-^'' can be represented as a sum of three matrices: 

.^(0) _ ^x{0) ^y{0) ^z(O) 
i] ij ~^ V i] ' 

2??^ = -(l/6)[<5..,,.-i + - '5..,,. (2 - <5,.,i - S,^.^m,^,,J,^,,^, (4.44) 

etc., each of them acting as a discrete Laplace operator along the respective coordinate. Thus the variables 
separate 

Fkt = Z^-.fcx X fiy,ky X f^^.k^, /^k = Afcx + i^ky + fj.k,, (4.45) 
and one obtains — 1 ^ ^k ^^'^ 

/i„,fc„ oc cos[(j„ - 1/2)A:J, fc„ = 7m„/iV, n^, = 0, 1, . . . , iV - 1 (4.46) 

with a = z,y, z. 
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Now, after substituting Eq. ( 4.421 ) in Eq. (4.41), one obtains for 6* ^ 1 



(fbc) 



A/" 4^ 1-; 



(4.47) 



[of. Eq. (3.33)]. It is interesting to note that VfI^^^ differs from W'^°'^^ of Eq. (3.21) only by the definition 



(pbc) 



the surfaces. In contrast to w'^'°'^\ the value of W^"'"'' exceeds the bulk value W , because the spins near 
the surface fluctuate more strongly at finite temperatures. Moreover, it can be shown that in the isotropic 
semi- infinite system the magnetisation profile is long ranged and decays as 1/^, where I is the distance from 



of the discrete wave vectors, Eqs. (3.17) and (4.46), which results from different boundary conditions at 



Kfbc) 



the surface. Integrating this profile over the particle's volume gives A 



(fbc) _ 
N — 



(fbc 
N 



W)/W cx ln(7V)/iV, 



i.e., the difference from the bulk value for the fbc system decreases more slowly than A^'^'^-' at large N . To 



obtain this result explicitly, one should bring Wf"^^ and I^^P''^) to the for ms that arc comparable with each 



other. To this end, we write W, 



(pbc) 



2N 



in the form 



(pbc) 



2N 



8M ^ 1 
k 



N 



,N-1 



(4.48) 



and complete the range of k in Eq. (4.47) to that of the formula above. This results in the exact relation 



(fbc 



N 



Wr 



(pbc 



2N 



167V 



E'- 



1 



(2)^ 



1 



r {i-x',>){i-\t>){b^\t>) 
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(4.49) 



where Aj^"* = {coskx +cosfcj,)/2, A^^"* = cosfc, and the additional terms arise because of the double counting 
on different faces, edges, and corners of the Brillouin zone as the completion is done. These terms can be 
loosely interpreted as contributions from faces, edges, and corners of the cube. One can see that in the 
large- iV limit the face contributions are the most important ones and they behave as \n{N)/N, whereas the 
edge contributions ^ 1/A^ and the corner contributions ~ 1/A^^. The result for A^'^'^-' has the form 



A 



(fbc) 
N 



w, 



(fbc 



N 



W _ 91n(1.177V) 



W 



2nNW 



iV> 1, 



(4.50) 



cf. Eq. ( 3.22 ). The next-to-log term in Eq. ( 4.50 ) is mo re di fficult to obtain analytically; the best way is to 
fit the result of the direct numerical calculation of Eq. (4.47) (see Fig. El). 



5 Numerical results 

Here we present our numerical results for the temperature and field dependence of the intrinsic magnetisation 
M and induced magnetisation m for small magnetic systems with simple cubic lattice and cubic shape, 
subject to free and periodic boundary conditions (fbc and pbc). As was mentioned in Sec. |[ the numerical 
method for solving the D oo model consists in obtaining the correlation function Sjj and magnetisation 



rrii from the linear equations (2.7), substituting them into the constraint equation (|2.6|), and solving the 



resulting nonlinear equation for the gap parameter Gi. On the first step we use a linear-equation solver 



based on the sparce-matrix storage appropriate for the Dyson matrix P, which is defined by Eq. (2.7), in our 
case of nearest-neighbour interactions. The nonlinear problem is tackled using a nonlinear-equation solver 
of the Newton-Raphson kind. For the present geometry, symmetry arguments have been used to reduce 
the number of unknowns and thereby to increase the computing performance. The relevant subset of sites 
corresponding to different values of Gi is obtained by taking 1/8 of the cube with a subsequent reduction by 
a factor of about 3! using permutations of the coordinates x, y, and z. This allows us to reduce the number 
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Fig. 1: Lattice sums Wn for the systems of the cubic form with free and periodic boundary conditions. 
W = 1.51639 is the bulk value for the sc lattice. 



of unknowns by a factor of about 8x3! = 48, this estimation becomes asymptotically exact for A/" ^ 1. 
For the square lattice the reduction factor is about 4 x 2! = 8. All computations have been performed on a 
Pentium III/450 MHz of 128 MB RAM. 

Fig. 1^ shows the temperature dependence of the intrinsic magnetisation A/, Eq. (2.11), and local mag- 
netisations Mi, Eq. (2.15), of the 14^ cubic system with free and periodic boundary conditions in zero field. 
For periodic boundary conditions, M exceeds the bulk magnetisation at all temperatures. In particular, at 
low temperatures this is in accord with the positive sign of the finite-size correction to the magnetisation, 
Eqs. (3.31) and (3.22). The magnetisation at the center of the cube with free boundary conditions is rather 
close to that for the model with pbc in the whole temperature range and converges with the latter at low 
temperatures. Local magnetisations at the center of the faces and edges and those at the corners decrease 
with temperature much faster than the magnetisation at the center. This is also true for the intrinsic mag- 
netisation M which is the average of the local magnetisation Mi over the volume of the system: For the 
relatively small size 14"^ the contribution from the boundaries to the average properties are still substantial. 
One can see that, in the temperature range below the bulk critical temperature, M is smaller than the bulk 
magnetisation. This means that the boundary effects suppressing M are stronger than the finite-size effects 
which l ead t o the increase of the latter. This is also seen from the low-temperature expresion for M given 
in Eq. (|l47|) . 

Magnetization profiles in the direction from the center of the cube to the center of a face are shown in 
Fig. H at different temperatures. One can see that perturbations due to the free boundaries extend deep 
into the system. This is a consequence of the Goldstone mode which renders the correlation length of an 
isotropic bulk magnet infinite below Tc- The latter effect is better seen from the magnified data at the low 
temperature (see the lower plot), where the MFA predicts, on the contrary, a fast approach to a constant 
magnetisation when moving away from the surfaces MC simulations of the classical Heisenberg model 
H, y also show long-range magnetisation profiles. The difference of the local magnetisations in the center 
and at a surface point reaches a maximun in the vicinity of the bulk Tc and goes to zero at low and high 
temperatures. 

Fig. ^ indicates that the critical indices for the magnetisation at the faces, edges, and corners are higher 
than the bulk critical index /3 — 1/2 for the present D = oo model. The critical index at the face pi is the 



mostly studied surface critical index (see, for a review, Refs. [g2, 23 ). The exact solution of Bray and Moore 
for the correlation functions at criticality in the D = oo model and application of the scaling arguments 
yield the value Pi — 1 (see Table II in Ref. P2|). Exact values of the edge and corner magnetisation indices, 
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Fig. 2: Temperature dependence of the intrinsic magnetization M, Eq. (2.11), and local magnetizations M^, 
Eq. {2.1L), of the 14"^ cubic system with free and periodic boundary conditions in zero field. 



(^2 and /33, seem to be unknown for D = oo. Cardy |^ used the first-order e-expansion to obtain /?2(a) for 
the edge with an arbitrary angle a. For a = 7r/2 and I? = cx3 in three dimensions the result for the edge 
critical exponent reads (32 = 13/8 + 0{e^) = 1.625 + O(e^). 

To estimate the magnetisation critical indices in our model we have performed a finite-size-scaling analysis 
(see, for a review, Ref. ||2^) assuming the scaling form M = N~^/'^Fm{tN^/'') and plotting in Fig. || the 
magnetisation times N^'^vs tN^/'^ . Here v = 1 is the critical index for the correlation length in the bulk 
and r = T jT^— 1, where = T^'^^ jW is the bulk Curie temperature. Our results for the systems with 
iV = 10 and = 14 merge into single "master curves" for /3i = 0.86, (52 = 1.33, and (i-i = 1.79, which have 
been obtained by fitting M cx iV-/?/'^ at T = i.e., 9 = Oc ^ 1/W. Note that our value 0.86 for the surface 
magnetisation critical index /3i is substantially lower than the value Pi — 1 following from scaling arguments. 
This disagreement is probably due to corrections to scaling which could be pronounced for our insufficiently 
large linear sizes TV = 10 and 14. A more efficient way for obtaining an accurate value of /3i is to perform 
a similar analysis for the semi-infinite model. The latter was considered analytically and numerically for 
T > Tc and ff = in Ref. jlj] . We also mention the Monte Carlo simulations of the Ising model |Q which 
yield Pi = 0.80, P2 = 1-28, and P3 = 1-77. 

The field dependence of M and m at fixed temperature, as obtained from the numerical solution of Eqs. 
(2.4) - (|2.6| ), for cubic and square systems is shown in Fig. ^. Naturally the numerical results for m confirm 
Eq. ( ^.13 ) which describes both the effect of orientation of the system's magnetisation by the field and the 
increase of M in field. On the contrary, using the zero-field value of M in Eq. ( ^.13| ) leads to a poor result 
for m shown by the dashed curve for the 10^ system in Fig. |[ Here we do see that there are two distinct field 
ranges separated by the characteristic field h* , which were discussed at the end of Sec. ||. The results of Fig. 
||are in qualitative agreement with those of Fig. 1 in Ref. jl8| showing the field dependence of x = dm/dh 
obtained from phenomenological arguments. We would like to stress, however, that Eq. (2.13) which is exact 
in the D = 00 model, has not yet been checked for systems with a finite number of spin components D and 
arbitrary size Af. The field dependence of the particle's magnetisation similar to that shown in Fig. || for 
the cubic system was experimentally obtained for ultrafine cobalt particles in recent Ref. as well as in 
a number of previous experiments. 

The curves for the square system in Fig. ^ illustrate the fact that in two dimensions thermal fluctuations 
are much stronger than in 3d, which leads to lower values of both M and m at the same temperature. 
The bulk magnetisation rrib in two dimensions vanishes at zero field and it thus goes below the intrinsic 
magnetisation M in the low-field region. 
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Fig. 3: Magnetization profiles in the direction from the center of the cube to the center of a face at different 
temperatures 9 = T/T^^^. The magnified lower plot shows the long-range character of the profile at low 
temperatures. 



6 Conclusion 

In this paper, we have studied finite-size and boundary effects in small ferromagnetic systems with free 
boundaries within the D ^ oo classical vector model. This model, while sacrificing about 5% in the overall 
accuracy when used as a substitute for the Heisenberg model in three dimensions, is exactly solvable and 
describes relevant physical features related to the Goldstone modes, such as the absence of ordering in two 
dimensions for isotropic models, long-range magnetisation profiles, etc. The D = oo model is clearly superior 
to the mean field approximation which ignores these important effects. 

An important result of our work is the exact relation in Eq. (2.13) beetween the intrinsic magnetisation 
M(H, T) and induced magnetisation m{H, T). It would be interesting to check whether this plausible relation 
holds for a more realistic Heisenberg model. In the latter case, it should be at least a good approximation. 

Whereas the finite-size effects, which occur in a pure form in systems with periodic boundary conditions, 
lead to the insrease of the magnetisation with respect to the bulk value, boundary effects in systems with 
free boundaries work in the opposite direction. Since the magnetisation reduction caused by the bound- 
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Fig. 4: Finite-size-scaling plot for the magnetization of cubic systems with linear size iV = 10 and 14. 



aries extends over the whole system (i.e., the magnetisation profiles are long ranged), the boundary effect 
overweighs the finite-size effect by a large logarithmic contribution (see Fig. |l|), and the resulting intrinsic 
magnetisation M is lower than the bulk magnetisation below the bulk critical temperature. The latter effect 
and the long-range magnetisation profiles are illustrated in Figs. ^ and 

We have also attempted to estimate the surface magnetisation critical indices for the faces, edges, and 
corners of the cube from the scaling analysis of the data for cubic systems with linear sizes = 10 and 
N = 14. It seems that our values for all these indices are too low and the data for larger sizes are needed. 
The required computer resources, however, rapidly increase with the system's linear size iV, since the theory 
operates with the correlation matrices of the size x A^'^. 

Field dependences of the particle's induced magnetisation m in Fig. |^ reflect both the alignment of the 
intrinsic magnetic moment M by the field and the increase of M in field. A combination of these two effects 
was observed in many experiments including the recent publication p8[ |. 

The present approach can be extended to magnetic systems of other shapes (spherical, cylindrical, etc.), 
more complicated lattice structures, and systems with a bulk and surface anisotropy. This will also allow 
us to compare with the results of the Monte Carlo simulations of a nanoparticle (of the maghemite type) 
studied in B. We are planning to study these issues in our future work. 
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